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The classical approximation may be applied to a number of problems in non-equilibrium field theory. The 
principles and limits of classical real-time lattice simulations are presented, with particular emphasis on the 
definition of particle numbers and energies and on applications to the earliest stages of heavy-ion collisions. 



1. INTRODUCTION 

One of the most important predictions of lat- 
tice QCD is the existence of a transition from 
ordinary hadronic matter, where quarks and glu- 
ons are confined, to a quark-gluon plasma at high 
temperatures. Current estimates are that this 
transition takes place at T c ~ 150 — 170 MeV, 
and is most likely a crossover at zero chemical 
potential. At higher chemical potential a first- 
order phase transition is predicted, ending in a 
tricritical point at ~ 100 - 500 MeV. 

These predictions are currently being put to the 
test at RHIC and other heavy-ion colliders. How- 
ever, interpreting the results of these experiments 
and comparing them to lattice QCD predictions is 
far from straightforward. One important reason 
for this is that the process to a large extent takes 
place out of thermal equilibrium, and equilibrium 
field theory methods such as lattice Monte Carlo 
are therefore not sufficient. 

One major, still unsolved puzzle of heavy- 
ion physics is whether the system ever actually 
reaches thermal equilibrium, and if so, what the 
equilibration time is. Many aspects of the colli- 
sion can be successfully described (eg, using hy- 
drodynamics P|) by assuming a very short equi- 
libration time of 1 fm/c, but it is far from un- 
derstood how this would come about. Clearly, a 
proper understanding of the thermalisation pro- 
cess is essential if we are to have a coherent and 
reliable description of the heavy-ion collision. 

At late stages of the collision process, as the 
system expands, the particles eventually decou- 
ple. This typically occurs in two phases: first, 
inelastic collisions cease, causing the ratios of dif- 



ferent particle species to be fixed (chemical freeze- 
out); later, also the mean free path for elastic 
collisions becomes too large and the momentum 
distribution of the particles deviates from that 
of thermal equilibrium (thermal or kinetic freeze- 
out). It is at this point that the final parti- 
cle yields of the collision are fixed. Thus, non- 
equilibrium dynamics is needed to understand 
also this aspect of the process. 

Another field where non-equilibrium field the- 
ory is needed, is early-universe physics. Examples 
of non-equilibrium processes which require a field- 
theoretical approach include (p)reheating, elec- 
troweak and QCD phase transitions in the early 
universe, and baryogenesis. 

Field theories out of equilibrium is a noto- 
riously difficult problem to study nonperturba- 
tively. A large number of approaches have been 
employed, including Hartree and large-N ap- 
proximations |2I3| . Dyson-Schwinger-related ap- 
proaches based on the 2PI effective action |4I5| . 
and kinetic theory 0. All of these have their 
strengths as well as drawbacks and limitations. 
Within its area of applicability, the classical ap- 
proximation has the advantage of being fully non- 
perturbative, easy and relatively inexpensive to 
implement numerically, and straightforwardly ap- 
plicable to gauge theories. The major drawback 
is obviously that quantum effects are not taken 
into account. 

2. THE CLASSICAL APPROXIMATION 

Classical statistical physics is quantum statis- 
tical physics in the limit of large occupation num- 
bers. It follows that the classical approximation 



can be used when the occupation numbers of the 
relevant or dominant modes of the system are 
large. In the context of heavy-ion collisions it 
is arguable that the multiplicity of soft gluons in 
the initial (pre-equilibrium) stages is very high, 
and the classical approximation may therefore be 
valid. 

The classical approximation has also been ap- 
plied to the chiral dynamics during freeze-out. In 
the linear sigma model and related models, the ef- 
fective potential of the chiral (a, 7?) field changes 
from symmetric to Mexican-hat type as the tem- 
perature drops, giving rise to an instability where 
the low-momentum modes increase exponentially. 
In this scenario, the resulting high occupation 
numbers justify the use of the classical approx- 
imation. 

A particular hazard with simulating classical 
dynamics on a lattice is connected with high- 
momentum lattice artefacts. These will in gen- 
eral interact with the soft modes which carry 
the interesting physics and for which the classi- 
cal approximation is in principle valid. If there is 
sufficient strength in the hard modes, they may 
equilibrate classically with the soft modes on a 
much shorter timescale than that of the interest- 
ing physics |7]. In that case, not only the hard 
modes, but the entire system will be dominated 
by classical lattice artefacts. 

To avoid this problem, it is important that the 
high-momentum modes should be, and remain, 
strongly suppressed. At early stages, this can 
be ensured by choosing appropriate initial con- 
ditions. 

The initial conditions are a crucial part of the 
simulation. They should reflect the salient fea- 
tures of the system at the outset. This is also 
the only place where information about the quan- 
tum nature of the real world enters into the sim- 
ulation. A quantum system may be represented 
as an ensemble of classical configurations initially 
distributed according to quantum statistics 
One example of this may be to choose the 2-point 
correlators of the fields and their canonical mo- 
menta to obey the Bose-Einstein distribution for 
free fields at some temperature T, after subtract- 
ing the quantum vacuum fluctuations. For scalar 
fields 4> with momentum fields it and mass m one 



would then have 

(0(fc,t = O)0(-fc,O)) = ^ — 7 i— , (i) 

(7r(fc,0)7r(-M)) =^k -^r— [ , (2) 

where u>k — \/fc 2 + m 2 . Such an initialisation 
also provides an exponential cutoff for the hard 
modes, which will help in avoiding the dangerous 
lattice artefacts. Following this, each configura- 
tion evolves independently according to the clas- 
sical hamiltonian equations of motion, and time- 
dependent correlators are computed as averages 
over the initial conditions. 

Given sufficient time, the system will even- 
tually thermalise classically, resulting in classi- 
cal equipartition = ^k/T and giving rise to 
Rayleigh-Jeans type divergences. The hope is 
that this will happen on much longer time scales 
than those under consideration in the simulation. 

3. DISTRIBUTION FUNCTIONS 

There is no unique definition of local particle 
numbers and energies for interacting fields out of 
equilibrium. Still, the system may exhibit effec- 
tive particle-like behaviour, which may be used 
to characterise the approach to thermal equilib- 
rium or to an equilibrium-like distribution. Given 
a definition of local particle numbers, these can 
also be used to give an effective description of the 
system in terms of kinetic theory. 

The effective particle numbers may be ex- 
tracted from the two-point field correlators, which 
in the free-field case (where the particle descrip- 
tion is appropriate and well-defined) contain all 
the information there is about the system. For 
example, for a homogeneous, free scalar field we 
have 

= -Ifn fc + ~) (3) 

(n{k,t)n(-k,t))=uj k (n k + ^J (4) 

where is the occupation number for the mode 
with momentum k and LUk is the associated en- 
ergy. For interacting fields, this may in turn be 



used as a definition of the instantaneous particle 
numbers and energies n and u> |9I8| : 



energies can be denned as 



n k (t)+- = VWMM-M)>MMM-M)> . 



u k (t) = >(MM-M)> 
V MMM-M)) 



(5) 
(6) 



In the classical approximation the 1/2 is left out. 

In a non-abelian gauge theory, the correlation 
functions will in general be gauge dependent, so 
the distribution functions will contain ambigui- 
ties due to the gauge choice. This ambiguity may 
be removed by constructing gauge invariant cor- 
relators using parallel transporters; however, this 
introduces path dependence. In particular, in a 
lattice regularisation there is in general no one 
preferred path between two points. 

Although the distribution functions are not 
unique, all physical observables extracted from 
them, such as masses, temperatures and chemi- 
cal potentials, should not depend on the defini- 
tion and in particular on the gauge. As long as 
these quantities are not well-defined on the other 
hand (such as when the system is very far from 
equilibrium and the quasiparticle picture does not 
apply), one may expect "masses" and "tempera- 
tures" to be definition-dependent. Thus, study- 
ing the gauge dependence (or path dependence) of 
distribution functions may serve the double pur- 
pose of monitoring the approach to equilibrium 
and verifying the validity of the approach used. 

One natural choice of gauge is the Coulomb 
gauge, which is a smooth gauge. In a system 
with spontaneously broken gauge symmetry (e.g., 
a Higgs system), the unitary gauge, where the (ef- 
fective) Higgs field has only one non-zero, real 
component, is another natural choice. Other 
gauges, such as maximal abelian gauge, axial 
gauges or random gauge, may also be considered. 

In the Coulomb gauge, the gauge potential A, 
(but not its conjugate momentum Ei) is purely 
transverse, and it can be shown that the trans- 
verse free correlators behave analogously to the 
scalar case ^UJ- Thus the particle numbers and 



n k = y/D^k)D^(k), u k = J£^Q. (7) 



Here D^,, D E are the transverse A- and i?-correla- 
tors respectively, constructed from the two-point 
functions ($£ (Jfc)$|(-jfe)) = 5 ah Cf?(k); $ = A, E 
according to 
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In figure ^ the gauge dependence of effective 
particle numbers is illustrated in the SU(2)-Higgs 
model In this case, the system was prepared 
in such a way that all the energy initially was 
in the Higgs field, while the gauge potential was 
initialised to zero. Since the angular modes of 
the Higgs fields are absorbed into the gauge fields 
in the unitary gauge, the inital occupation num- 
bers are very different. However, already after 
t = 10 rnTj the two distributions appear almost 
identical. However, while the particle numbers in 
the Coulomb gauge change very little from here 
on, in the unitary gauge they continue to fluc- 
tuate and it is only from f « 40 m ~H orL that 
one with some confidence can claim the numbers 
are gauge independent. This agrees roughly with 
the point where the dispersion relation in the uni- 
tary gauge begins to show stable, particle- like be- 
haviour in this particular model. 

In a non-homogeneous system, and in general 
in a kinetic-theory description, it is appropriate to 
think in terms of local particle numbers n(x, k, t). 
These may be related to the Wigner functions 
constructed from gauge invariant two-point func- 
tions 9 , or more generally to the two-point func- 
tions fourier transformed on a region R(x) cen- 
tred on x. In the case of a scalar theory we may 
have (suppressing the common t-coordinate for 



brevity) 
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Figure 1. Gauge dependence of gauge particle 
numbers in the SU(2)-Higgs model ^U] a t differ- 
ent times. The filled symbols are particle num- 
bers obtained in the Coulomb gauge, while the 
open symbols are for the transverse modes in the 
unitary gauge. 



il R JR{x) 

(10) 

C^(x,k) = -1- / d 3 zd 3 ye~^-^(n(y)ir(z)) . 

* l R JR{3) 

(11) 

Here is the volume of the region R(x). This 
coarse-graining creates an intrinsic unsharpness 
in the momentum k and position x of the quasi- 
particles, given by the size (and shape) of the re- 
gion R. 

If the system under consideration is homoge- 
neous, we may improve statistics by performing 
an average over all space. This can be shown 
to be equivalent to local averaging in momentum 
space, with a weight function w depending on the 
size and shape of R: 

C{k,t) = i Jd 3 xC(x, k,t) 

= — Jd 3 xC R{S) {k,t) (12) 

= ^2 W (P- k)C v (p,t) , 
p 

where Cy denotes the correlation function evalu- 
ated on the total volume V, and the sum is over 
the discrete momenta available on this volume. In 
practice, it is simpler to work backwards, choos- 
ing a simple form of momentum-space averaging 
which may correspond to rather complicated spa- 
tial regions. For instance, binning in the absolute 
value of the momentum, 

w (f-k)ocQ(\k\-\0\+^)-QQk\-\p\~), (13) 

corresponds to spherical shells with thickness ap- 
proximately 1/A in position space. 

4. PURE YANG MILLS 

At the earliest stages of heavy-ion collisions, 
the gluon density is expected to be so high that 
the classical approximation can be justified. The 
same approximation also justifies ignoring the 



back-reaction of the quarks, since their number 
density will be much lower; leaving us with clas- 
sical Yang-Mills equations of motion, which may 
be solved numerically on a lattice. 

The lattice equations of motion in the temporal 
gauge (Aq — 0) read 

dtE?(x)=Dftv\it b U ji (x)] (14) 
where 

E?(x) = F-(x) = tv^U^Ujix + t)] (15) 

is the canonical momentum to Af. Here, dt de- 
notes the backward lattice derivative, while Dj 
is the backward covariant lattice derivative. The 
equations of motion for Aq constitute the Gauss 
constraint, 

D?E$ = 0, (16) 

which must be satisfied by the initial conditions 
but is conserved by the equations of motion. 

The initial gluon fields should be related to the 
gluon distributions of the two colliding nuclei: in 
principle they should just be the superposition of 
two Lorentz-boosted nuclear gluon distributions. 
Simulations have been carried out over a num- 
ber of years by Krasnitz, Nara and Venugopalan 
|11I12| (see also ^1) using the "colour glass con- 
densate" model of the nuclear wave function to 
provide the initial conditions. In these studies, 
the numerical work has been simplified by con- 
sidering only the mid-rapidity region where the 
physics is assumed to be boost-invariant. This 
reduces the system to effectively 2+1 dimensions. 
With these assumptions, the authors have been 
able to provide an estimate of the initial energy 
density and gluon distribution which may be used 
as input into hydrodynamic or kinetic calcula- 
tions. 

An alternative approach would be to determine 
the nuclear gluon field from e.g. a bag model, 
give this a Lorentz boost, and perform a 2+1+1- 
dimensional simulation with the longitudinal lat- 
tice spacing a z = aj/7, where a± is the lattice 
spacing in the transverse (cc, y) direction. Work 
is underway to implement this. 



5. SUMMARY 

The classical approximation may be applied 
to a range of problems in non-equilibrium field 
theory where occupation numbers are high, such 
as the earliest stages of heavy-ion collisions. It 
has the advantage of being non-perturbative and 
computationally relatively inexpensive. Effective 
particle numbers may be defined out of equilib- 
rium in a self-consistent manner, and their gauge 
dependence (or that of derived quantities such as 
masses and temperatures) can be used as a check 
on the validity of the quasi-particle picture. 
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